Matlab code 3.3: Matlab file “Figure 3-12.m”
%--------------------------------------------------------------------
% This
code can be used to generate Figure 3.12
%--------------------------------------------------------------------
clear
all
close
all
%% the
theoretical curve of the precession target micro-Doppler characteristic(LFM
radar)
c = 3e8;
j =
sqrt(-1);
fc =
20.5e9; % carrier frequency of transmitted signal
tp =
50e-6; % pulse duration
B = 3e9;
% bandwidth
prf =
1000; % pulse repetition frequency
pri =
1/prf; % pulse repetition interval
t = 1.5;
% duration of echo signal
cord =
1000*[-100 -300 500]; % coordinates of local coordinate system's origin in the
radar coordinate system
colo =
[0 0 0];
v = 0; %
translational velocity of target
ws =
pi*[0 0 6]; % angular velocity of spinning in the local coordinate system
omegas =
6*pi;
Ts =
0.33; % spinning period
wc =
[5.4553 0 30.9386]; % angular velocity of coning in the local coordinate system
omegac =
10*pi; % angular frequency
Tc =
0.2; % coning period
p1 = [0
0 2]; % scatterer on the apex of cone
p2 =
[-0.6 0 -1]; % scatterer on the bottom of cone
p3 =
[-0.6 0 1]; % scatterer on the bottom of cone
tgt =
[p1;p2;p3];
ae =
pi*[1/3 1/4 1/5]; % Euler angle
ri1 =
[cos(ae(3)) sin(ae(3)) 0;-sin(ae(3)) cos(ae(3)) 0;0 0 1];
ri2 = [1
0 0;0 cos(ae(2)) sin(ae(2));0 -sin(ae(2)) cos(ae(2))];
ri3 =
[cos(ae(1)) sin(ae(1)) 0;-sin(ae(1)) cos(ae(1)) 0;0 0 1];
ri =
ri1*ri2*ri3; % initial rotation matrix
ws =
ri*ws'; % angular velocity of spinning in the reference coordinate system
wse =
ws/omegas; % unit vector of spinning
wsr = [0
-wse(3) wse(2);wse(3) 0 -wse(1);-wse(2) wse(1) 0]; % skew symmetric matrix of
spinning
wc =
ri*wc'; % angular velocity of coning in the reference coordinate system
wce =
wc/omegac; % unit vector of coning
wcr = [0
-wce(3) wce(2);wce(3) 0 -wce(1);-wce(2) wce(1) 0]; % skew symmetric matrix of
coning
r0 = ri*tgt;
% scatterers in the reference coordinate system
R0 =
colo-cord;
n =
R0/sqrt(sum(R0.^2)); % unit vector of Line-of-Sight(LOS)
dt =
0:pri:t-pri; % slow time sampling interval
m =
length(dt);
Rtm =
zeros(m,length(tgt));
for i =
1:m
Rc = eye(3)+wcr*sin(omegac*dt(i))+wcr^2*(1-cos(omegac*dt(i)));
Rs =
eye(3)+wsr*sin(omegas*dt(i))+wsr^2*(1-cos(omegas*dt(i)));
for k = 1:length(tgt)
if k == 1
Rtm(i,k) = r0(k,:)*ri'*Rc'*n';
else
Rtm(i,k) = r0(k,:)*ri'*(Rc*Rs)'*n';
end
end
end
figure(1)
plot(dt,Rtm)
xlabel('Slow
time (s)')
ylabel('Range
(m)')
axis([0,1.5,-2.5,2.5])
%% the
simulation result of the precession target micro-Doppler characteristic(LFM
radar)
rmin =
sqrt(sum((cord-colo).^2)); % the minmum detectable distance
rmax =
sqrt(sum((cord-colo).^2))+sqrt(sum((tgt(1,:)-colo).^2)); % the maximum
detectable distance
ts =
2*rmin/c; % the starting point of the time sampling
te =
2*rmax/c+tp; % the ending point of the time sampling
mu =
B/tp; % frequency modulation rate
fs =
B/100;
kdt =
1/fs; % fast time sampling interval
nt =
2*ceil((te-ts)/(2*kdt)); % number of fast time sampling
df =
fs/nt; % frequency sampling interval
Rr =
zeros(m,length(tgt));
for i =
1:m
for k = 1:length(tgt)
Rr(i,k) = sqrt(sum(R0.^2))+v*dt(i)+Rtm(i,k);
end
end
s =
zeros(nt,length(dt)); % echo signal
tk =
ts+kdt*(0:nt-1); % fast time
for n =
1:length(tgt)
td =
tk(:)*ones(1,m)-2*ones(nt,1)*(Rr(:,n))'/c;
s = s+exp(j*2*pi*fc*td+j*pi*mu*td.^2);
end
Rr0 =
sqrt(sum(R0.^2));
td0 =
tk(:)*ones(1,m)-2*Rr0/c;
s0 =
exp(j*2*pi*fc*td0+j*pi*mu*td0.^2); % reference signal
sdt =
s.*conj(s0); % dechirp
sdf =
fftshift(fft(sdt),1);
x = dt;
y =
-((-nt/2:nt/2-1)*df)*c/(2*mu);
figure(2)
contour(x,y,abs(sdf),15)
xlabel('Slow
time (s)')
ylabel('Range
(m)')
axis([0,1.5,-2.5,2.5])